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Abstract 

We present a practical application of parallel symbolic computation in General Relativity: the 
calculation of curvature invariants for large dimension. We discuss the structure of the calculations, 
an implementation of the technique and scaling of the computation with spacetime dimension for 
various invariants. 
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Parallel symbolic computation has received increasing attention since the early 1990's, 
and a number of applications in pure mathematics and databases have been ex- 

plored. Little attention, however, has been paid to possible applications in General Relativ- 
ity. This note describes such an application. 

Curvature invariants are scalar products of Riemann, Ricci or Weyl Tensors or their 
covariant derivatives. The canonical example is the Kretschmann Invariant 

-K J^abcd 

Here Einstein Summation (also known as index contraction) is performed on all repeated 
indices. The Kretschmann Invariant is the simplest invariant product involving the Riemann 
Curvature Tensor. It is used most often to identify essential singularities in a spacetime 
geometry. 

A variety of other curvature invariants have been considered in the literature. For example 

II — n nagch;ef-n, rtibjd;kl 

(where the semi-colon denotes covariant differentiation) has been used P] to analyze type N 
spacetimes. The same authors jsj used 

h = R"'^'^'^RaefgR^^ bhR^^ cd (3) 

to locate asymptotically flat regions of the spinning C metric. These calculations all share 
two common characteristics: 

• they are sums of products of tensor components whose indices can all be enumerated 
before the calculation begins, and 

• the final results are relatively simple sums of multinomials. 
For example, I2 for the spinning C metric is P| 

-144™',^^. (4) 



Since the set of indices occurring in the sum is enumerable before computation begins, 
the sum can easily be partitioned and carried out on multiple processors. Because the metric 
is typically a function of a small number of variables, the number of different multinomials 
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occurring in the sum is correspondingly small relative to the number of products which must 
be summed. This means that partial sums can be accumulated on each processor, and when 
all of the partial sums are complete, final summation is in general a short process. For the 
same reason, simplification of the final sum is relatively straightforward. 

This application of the Divide and Conquer Principle yilends itself well to loosely coupled 
multiprocessor configurations such as the Beowulf Cluster 7| because the close similarities 
between the partial sums lead to an essentially linear scaling of execution time with cluster 
size. It is important to note that as the number of processors increases, the sizes of the 
partial sums may increase as terms are accumulated in different processors which will cancel 
in the final sum. This increase is in general small, but is dependent on the specific metric 
and invariant studied. Since the sizes of the various partial sums are not known a prioi, a 
slight overall decrease in execution time can be achieved by partitioning the workload into 
many more parcels than processors. 

In D dimensions, each summation index can take D values, so that an invariant with n 
summed indices requires the summation of at most products. This number can often be 
reduced by consideration of the symmetries of the Riemann Curvature Tensor: 

Rabcd Rbacd 
Rabdc 

= Rcdab (5) 

The first two of these symmetries implies that pairs of antisymmetric indices which are 
summed together (for instance, c and d in I2) contribute a factor of ^^^^^^"^ to the number of 
products. The last symmetry is only useful for specific invariants such as the Kretschmann 
Invariant, where it reduces the number of products by the factor 

(6) 

D{D-l) + 2 ^ ' 

Summation over pairs of symmetric indices (such as e and f or k and / in /i) contribute a 
factor of ^^^^^^^ to the number of products. Sums over lone indices contribute a factor of 
D. It must also be remembered that for any given metric, not all of the Riemann Tensor 
components are nonzero. Hence the computation of Ji involves the sum of at most 

4 ^ ^ 
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products quartic in second derivatives of the Riemann Tensor, each of which may contain 
tens of terms. In four dimensions, this is 6,553,600 products and perhaps tens or hundreds 
of biUions of terms, depending on the complexity of the metric. In comparison, the worst 
case scenario for the computation of I2 requires the sum of 

^ (8) 

products quartic in the Riemann Tensor. 

The algorithm for index enumeration is straightforward if a bit messy. Assume that 
tensor -array is an array of tensors occuring in the invariant and that all necessary index 
raising has already been accomplished. For instance, in I2, the elements of tensor -array 
would be 

TDabcd D D06 cnA T^oh 

-n- ) ^abcdi -K cd a-HCl K cd- 

The tensor data structure must contain, in addition to a sparse array containing the nonzero 
components for that tensor, an array indicating which pairs of indices are antisymmetric, 
and an array indices which indicates which indices are summed together. For /2, indices 
for each tensor would be 

{1, 2, 3, 4}, {1, 5, 6, 7}, {5, 6, 2, 8} and {7, 8, 3, 4}, 

indicating that, for instance, the second index of the first tensor should be summed over 
with the third index of the third tensor, etc. 

The following pseudocode details the enumeration of the indices of the component prod- 
ucts to be summed with respect to antisymmetric index pairs as well as the elimination of 
products involving zero components. We begin by figuring out which pairs of sums can be 
abbreviated due to symmetry considerations: 

multiplier — 1 

for each pair of tensors A and B { 

for each pair of indices{I, I + 1) in A and (J, J + 1) in B { 

if A.indices{I) == B.indices{J) and A.indices{I + 1) == B.indices{J + 1) { 
if both pairs are antisymmetric { 

remember to abbreviate summation on this pair of indices 
multiplier — multiplier * 2 } } } } 
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Summation abbreviation is accomplished by ignoring index pairs in which the second index 
of each pair is less than the first. Each product will be multiplied by multiplier to correct 
for the abbreviation. We chose not to implement summation abbreviation for symmetric 
pairs of indices, or for larger sets of antisymmetric indices, for reasons of simplicity. We 
also chose to ignore the third symmetry property of the Riemann Tensor, because it applies 
relatively rarely in the types of invariants which might benefit from parallel processing. 

We can now enumerate the indices which will contribute to the invariant. The array 
sumJndices contains the set of indices under consideration for a given product. We assume 
that the range of indices is from to D — 1, and that as a given index increments to D, it 
is set to zero and the next index is incremented by 1. In the following we will refer to this 
procedure as "cycling the indices". The resulting sumJndex -array is the complete list of 
all products which contribute to the invariant: 

set sumJndices to all zeroes 
product-Count = 
do { 

process JheseJndices — TRU E 
for each index pair {I, J} { 

if abbreviating summation on this pair { 

process JheseJndices = process dhescAndices AND [J > I) } } 
if process -these Jndices { 

if all tensor component corresponding to these indices are nonzero { 
sum Jndex -array (product -Count) — sumJndices 
product-Count — product-Count + 1 } } 
cycle sumJndices } 
until sumJndices cycles back to all zeroes 

The parallelization of curvature invariant computation can then be summarized as follows, 
assuming n processors with m parcels per processor: 

1. compute the components of the Riemann Tensor and any derivatives which occur in 
the invariant; 

2. enumerate the tensor indices occurring in the sum as described above; 
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3. partition sum Judex -array into m* n parcels; 

4. distribute new parcels to processors as they become idle, with each processor accumu- 
lating a partial sum; 

5. retrieve the n partial sums and simplify the full sum. 

This algorithm has been implemented by the author in a prototype code set called PTAH 

n 

(Parallel Tensor Algebra Hybrid) [9], which is a hybrid of C++ and Mathematica. The 
results below were obtained using that code, and were then cross-checked using Mathemat- 
ica alone on a single processor. The enumeration pseudocode above is derived from the 
corresponding code in PTAH. 

To examine the scaling properties of this technique, we require a related class of spacetime 
metrics in diverse dimensions. One such class is that of the Kerr Metrics with single rotation 
1^. Let 



parameter 



2 2,2 2 a 

p = r + a cos V 



fyD 5 

= TT—4: — — (9) 

(H + a^) - 

where fi and a are proportional to the mass and angular momentum respectively. Then the 
metric in D dimensions is 

ds^ = cos^ edQ^ + "^Ddr^ + p^de^+ 

((r^ + a^) sin^ 9 + A^ja^ sin^ e)d(f)^+ 

2 Ana sin^ edcpdt + (A^ - l)df (10) 

where dfl'^ is the standard metric on S*^"^. We note that the nonzero Riemann Tensor 
components for this metric in D = 4 have between 2 and 18 terms. In D = 6 and 8, 
the Riemann Tensor components have as many as 25 terms. For a general metric in D 
dimensions, there are 



(11) 



12 

independent curvature components. For the Kerr metric in D = 4, 13 of the possible 21 
independent components are nonzero. In contrast, for D = 11, the Kerr curvature has only 
68 nonzero independent components of the possible 1210. 



The Kerr Metrics depend explicitly on D — 1 variables: r, 9, a, fi and the D — 5 "polar" 
angles parametrizing the S^~^. We summarize the computations of three invariants: 

T T^abcd r> 

la — rC rCabcd 

T T^abcd -Qef r> 

Ih — n n ab^cdef 

T Tjabcd;e TD /"lO^ 

Ic — rt rtabcd;e [i-^) 

for D = 4 to D = 11 in the following table Q: 



D 
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1160 
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290 


391 


2749 
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21 
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373 


495 


3982 
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21 
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477 


628 


5703 
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21 
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606 


798 


8109 
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21 


10 


764 


1014 


11390 


5 


7 


21 


11 


955 


1286 


15769 


5 


7 


21 



Here P{I) is the number of products of tensor components in the sum and T(J) is the number 
of terms in the final result. The P{I) includes those operations necessary to raise tensor 
components before index contraction is possible, but excludes all unnecessary products due 
to zero components and symmetry considerations. 

Even with the rise in the number of symbolic multiplications with spacetime dimension 
and the substantial variation between different invariants, we see that for this class of ex- 
amples the final sums are quite small. This indicates that this technique can be efficiently 
parallelized under a loosely coupled MIMD architecture. We have also investigated a set of 
supergravity brane metrics in 10 and 11 dimensions which produced even smaller final 
sums (in accordance with the relatively simpler metric structures involved). 

While they illustrate the salient features of the computational technique, the Kerr invari- 
ants discussed here are not sufficiently complicated to warrant parallelization on a modern 
Beowulf Cluster. A more challenging computation is the Euler Class, which in 10 dimensions 
is quintic in the Riemann Tensor and requires (worst case) the sum of over 6 billion products. 
For the Kerr metric, the result has but 4 terms. For sufficiently complicated metrics and 
large D, invariants like Ji and I2 may not be practically computable in any other way but in 



parallel, but because of the small number of variables parametrizing the metric, the results 
are expected to be similar to those described here. 

It is also possible to use this technique to perform index contractions in tensor products 
which are not scalars. In this case the contracted indices are enumerated once for each 
set of values of the free (uncontracted) indices. The technique is even more effective in 
spinor contractions. Since the number of spinor components in D dimensions is 2^^"'"'^/^, 
the number of products occurring in each sum grows much more quickly with dimension 
than in the tensor case. 
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